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ABSTRACT 

Dwarf spheroidal (dSph) galaxies are prime targets for present and future 7 -ray telescopes 
hunting for indirect signals of particle dark matter. The interpretation of the data requires 
careful assessment of their dark matter content in order to derive robust constraints on can¬ 
didate relic particles. Here, we use an optimised spherical Jeans analysis to reconstruct the 
‘astrophysical factor’ for both annihilating and decaying dark matter in 21 known dSphs. Im¬ 
provements with respect to previous works are: (i) the use of more flexible luminosity and 
anisotropy prohles to minimise biases, (ii) the use of weak priors tailored on extensive sets 
of contamination-free mock data to improve the conhdence intervals, (iii) systematic cross¬ 
checks of binned and unbinned analyses on mock and real data, and (iv) the use of mock data 
including stellar contamination to test the impact on reconstructed signals. Our analysis pro¬ 
vides updated values for the dark matter content of 8 ‘classical’ and 13 ‘ultrafaint’ dSphs, with 
the quoted uncertainties directly linked to the sample size; the more flexible parametrisation 
we use results in changes compared to previous calculations. This translates into our ranking 
of potentially-brightest and most robust targets—viz., Ursa Minor, Draco, Sculptor—, and of 
the more promising, but uncertain targets—viz., Ursa Major 2, Coma—for annihilating dark 
matter. Our analysis of Segue 1 is extremely sensitive to whether we include or exclude a 
few marginal member stars, making this target one of the most uncertain. Our analysis illus¬ 
trates challenges that will need to be addressed when inferring the dark matter content of new 
‘ultrafaint’ satellites that are beginning to be discovered in southern sky surveys. 

Key words: astroparticle physics — (cosmology:) dark matter — Galaxy: kinematics and 
dynamics — 7 -rays: general — methods: miscellaneous 


1 INTRODUCTION 

Constraining particle candidates for dark matter (DM) through in¬ 
direct searches for their annihilations or decays has a long history 
(Gunn et al. 1978; Stecker 1978; for a recent review see Bergstrom 
2012). Several astrophysical messengers (radio/infrared; X-rays; 
7 -rays and neutrinos) and targets (e.g.. Galactic centre, dwarf 
spheroidal galaxies, galaxy clusters) have been used for such stud- 
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ies (e.g. Boyarsky, lakubovskyi & Ruchayskiy 2012; Cirelli 2012; 
Hooper 2012) with inconclusive results so far. 

Nearby dwarf spheroidal (dSph) galaxies are particularly 
promising targets because of their proximity and low intrinsic 
astrophysical backgrounds (Lake 1990; Evans, Ferrer & Sarkar 
2004). Indeed, several upper limits have been set on the thermally- 
averaged self-annihilation cross section {ov) of constituent dark 
matter particles using 7 -ray observations of dSphs. The H.E.S.S. 
collaboration (Abramowski et al. 2014) in their most recent anal¬ 
ysis which combines 5 dSph galaxies require {av) < 3.9 x 
10“^"^ cm® s“® (95% CL) for ©(TeV) mass particles. By ob- 
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serving Segue I for 160 hours, the MAGIC collaboration infers 
{(jv) < 1.2 X cm® s“® for ~ 500 GeV mass DM anni¬ 

hilating into the t'^t~ channel (Paiano et al. 2011; Aleksi et al. 

2014) .^ These values are to date the strongest limits set by ground- 
based imaging air Cherenkov telescopes. At lower energy, obser¬ 
vations of several dSphs by the Fermi-LAT satellite (Geringer- 
Sameth, Koushiappas & Walker 2014; Fermi-LAT Collaboration 
2014, 2015) provide the most stringent constraint to date with 
(cru) < 3 X 10“®® cm® s“® (95% CL) for dark matter particles 
with a mass below 100 GeV annihilating into the bb channel. As this 
is of the order of the thermally-averaged annihilation cross section 
required for a weakly interacting massive particle to make up the 
dark matter, it has become imperative to scrutinise these constraints 
closely. Intriguingly, Geringer-Sameth et al. (2015) report evidence 
for 7 -ray emission from the recently discovered Milky Way satel¬ 
lite Reticulum II (Koposov et al. 2015; DES Collaboration et al. 

2015) consistent with annihilating DM of mass 0{1Q — 100) GeV 
scale. The subsequent analysis of Hooper & Linden (2015) has con¬ 
firmed the apparent 7 -ray excess while Fermi-LAT Collaboration 
et al. (2015) claim no statistically significant detection. An under¬ 
standing of Reticulum IPs DM content will be crucial for deter¬ 
mining whether this interpretation of the 7 -ray signal is compatible 
with non-detections in other dSphs. 

In the X-ray band, an unidentified excess at 3.55 keV detected 
by stacking® XMM-Newton galaxy cluster spectra (Bulbul et al. 
2014a), and by using deep observations of M31 and the Perseus 
cluster (Boyarsky et al. 2014b) has been interpreted as a possible 
signature of decaying dark matter. However, this interpretation has 
been disputed (Park, Kong & Park 2014; Jeltema & Profumo 2014; 
Boyarsky et al. 2014a; Bulbul et al. 2014b). Notably, Malyshev, 
Neronov & Eckert (2014) find no such excess when performing 
a similar analysis on dSphs, possibly ruling out the decaying DM 
interpretation.® 

As mentioned, dSph galaxies have been highlighted as targets 
because of their high mass-to-luminosity ratio and absence of astro- 
physical emission at X-ray and 7 -ray energies (Lake 1990; Evans, 
Ferrer & Sarkar 2004). However to robustly use these objects to 
constrain particle DM candidates requires reliable estimates of the 
astrophysical factor for a potential signal (the so-called J- and 
79-factors, for annihilation and decay respectively). This requires 
careful and optimal evaluation of the DM distribution in these ob¬ 
jects. While cosmological priors or assumed DM profiles are of¬ 
ten used to this end (Sanchez-Conde et al. 2007; Fieri, Lattanzi & 
Silk 2009; Martinez 2013), data-driven approaches have been de¬ 
veloped, yielding more reliable results (Strigari et al. 2007; Essig, 
Sehgal & Strigari 2009; Geringer-Sameth, Koushiappas & Walker 
2015). 

We used such a data-driven approach earlier (Charbonnier 
et al. 2011) to obtain conservative estimates of the astrophysical J- 
factors of the 8 ‘classical’ dSphs — those with the best-measured 
stellar kinematics to date. Here, we revisit this question using the 
optimised kinematic data analysis developed by Bonnivard et al. 
(2015) to obtain reliable values for both J- and 79-factors for these 
8 ‘classical’ dSphs. Additionally, we include 13 ‘ultrafaint’ dSph 

^ The annihilation channel determines the energy spectrum of the 7 -rays, 
hence the instrumental response. 

^ Stacking signals from different galaxy clusters is an effective approach 
for decaying dark matter (Combet et al. 2012) but less so for an annihilating 
particle (Nezri et al. 2012; Maurin et al. 2012). 

® This non-detection has motivated new models for DM attempting to rec¬ 
oncile both claims, e.g. Cline & Frey (2014). 


galaxies in our analysis. Despite their sparse and often very uncer¬ 
tain kinematic data, the close proximity of some of these dSphs 
(Segue I is located at only ~ 23 kpc) has focused attention on 
them, and the most stringent limits set on {av) have used such 
objects. This makes it vital to obtain a proper measure of both sta¬ 
tistical and systematic errors for the astrophysical factors J and 79 
for these dSphs (a recent case in point is the tiny Hercules galaxy 
whose mass estimate dropped by a factor ~ 3 once contaminating 
foreground stars were properly weeded out; Aden et al. 2009). Es¬ 
pecially for pointed observations where many hours of integration 
may be lavished on a single target, it is vital to know not only which 
dSphs are the most promising targets, but also which are the most 
robust targets. This is our key goal: to obtain unbiased J and 79 
measures with robust uncertainties. 

The nearby dSphs with good data favour dark matter cores 
over the cusps predicted by pure cold dark matter cosmological 
simulations (e.g. Kleyna et al. 2003; Goerdt et al. 2006; Battaglia 
et al. 2008; Walker & Penarrubia 2011; Cole et al. 2012; Agnello 
& Evans 2012; Amorisco, Agnello & Evans 2013), though this has 
been disputed at least for some dSphs (e.g. Breddels & Helmi 2014; 
Richardson & Fairbairn 2014; Strigari, Frenk & White 2014). Such 
cores can result from stellar feedback from galaxy formation. This 
can be surprisingly effective even for low stellar mass systems (e.g. 
Navarro, Eke & Erenk 1996; Read & Gilmore 2005; Pontzen & 
Govemato 2012; Teyssier et al. 2013; Madau, Shen & Governato 
2014; Onorbe et al. 2015). Given such theoretical and observational 
uncertainties, in this paper we allow for considerable freedom in the 
central dark matter density and slope. 

This paper is organised as follows. In Sect. 2, we present the 
Jeans analysis and the input parametrisations for the DM, light, 
and anisotropy profiles. We discuss the several likelihood functions 
tested to match the data and introduce the astrophysical factors J 
and 79 calculated from the reconstructed DM profiles. Sect. 3 de¬ 
scribes the surface brightness and kinematic data used in the anal¬ 
ysis. In Sect. 4, we present the setup used in the Markov Chain 
Monte Carlo (MCMC) analysis, in particular our choice for the pri¬ 
ors of the free parameters. We also discuss the choice to be made 
for the DM halo size and the possible contamination of the data. 
Results are presented in Sect. 5, focusing on the overall contrast 
of the signals from dSphs with respect to that from the Galactic 
DM halo, and on the ranking of J- and 79-factors, in comparison to 
previous works. Conclusions are presented in Sect. 6 , while appen¬ 
dices A, B, and C provide further details regarding specific points 
addressed in this work. 


2 JEANS ANALYSIS, LIKELIHOOD FUNCTIONS AND 
ASTROPHYSICAL J- AND 79-FACTORS 

2.1 Jeans analysis 

The mass density profile of the dSphs is the key input for determin¬ 
ing their J- and 79-factors (Section 2.3). Different techniques have 
been developed in order to infer the mass profiles from stellar kine¬ 
matic data, such as distribution function modelling, Schwarzschild 
and ‘Made-To-Measure’ methods, as well as Jeans analysis (see re¬ 
cent reviews by Battaglia, Helmi & Breddels 2013; Strigari 2013; 
Walker 2013). In this work we focus on the latter, using para¬ 
metric functions as ingredients of the spherically symmetric Jeans 
equation. This technique has already been widely applied to dSphs 
(Strigari et al. 2007; Essig et al. 2010; Charbonnier et al. 2011; 
Geringer-Sameth, Koushiappas & Walker 2015). Here, we apply 


Dark matter annihilation and decay in dSphs 3 


the findings of Bonnivard et al. (2015), where an optimised strategy 
was proposed to mitigate possible biases introduced by the Jeans 
modelling. 


2.1.1 Spherical Jeans equation 


Dwarf spheroidal galaxies are considered as collisionless sys¬ 
tems described by their phase-space distribution function, which 
obeys the collisionless Boltzmann equation. Assuming steady- 
state, spherical symmetry and negligible rotational support, the 
second-order Jeans equation is obtained by integrating moments of 
the phase-space distribution function (Binney & Tremaine 2008): 

1 d - = ( 1 ) 

u dr r 


where v(r), v1{r), and /3ani(f) = 1 — are the stellar number 

density, velocity dispersion, and velocity anisotropy, respectively. 
Neglecting the (< 1%) contribution of the stellar component, the 
enclosed mass at radius r can be written as 

nr 

M(r) = 47 r/ pDM(s)s^ds, ( 2 ) 

^0 

where pDM(r) is the DM mass density profile. The solution to the 
Jeans equation relates M(r) to v(r)v‘l{r). However, the internal 
proper motions of stars in dSphs are not resolved, and only line-of- 
sight projected observables can be used: 


]{R) = 




a 


1 /3ani(r) ^2 


v{r) iip(r) r 
s/r'^ — 


dr, (3) 


with R the projected radius, o-p{R) the projected stellar velocity 
dispersion, and E(J?) the projected light profile (or surface bright¬ 
ness) given by 


E(i?) = 2 


Jr V 


(r) r dr 
s/r^ — R'^ 


(4) 


Note that the velocity anisotropy /3ani(r) cannot be measured di¬ 
rectly, in contrast to crp(J?) and E(J?). In our approach, parametric 
models for Pam{r) and pDM(r) are assumed in order to compute 
(jp (R) via equation (3). We can then determine the parameters that 
reproduce best the measured velocity dispersion crobs(^?)- 


2.1.2 Choice of parametric functions 


DM density profile: Following Charbonnier et al. (2011), we 
do not use a strong cosmological prior (e.g. assume the profile to 
be cuspy), as this will bias the derived astrophysical factors. In¬ 
stead, we fit the model parameters to data. We adopt the Einasto 
parametrisation of the DM density profile (Merritt et al. 2006): 


Einasto 

Pdm 


(r) 


= p- 2 exp 


2 

a 



(5) 


where the three free parameters are the logarithmic slope a, the 
scale radius r _2 and the normalisation p- 2 . Bonnivard et al. 
(2015) find that the choice of parametrisation — Zhao-Hernquist 
or Einasto — has negligible impact on the calculated J- or D- fac¬ 
tors and their uncertainties. With fewer free parameters, the Einasto 
parametrisation is more optimal in terms of computational time. 

Velocity anisotropy profile: We use the Baes & van Hese 
(2007) parametrisation to describe the velocity anisotropy profile: 


where the four free parameters are the central anisotropy /3o, the 
anisotropy at large radii /Joo, and the sharpness of the transition rj 
at the scale radius ra. This parametrisation was found to mitigate 
some of the biases arising in the Jeans analysis when using less 
flexible anisotropy functions with fewer free parameters (e.g., con¬ 
stant, Osipkov-Merrit — see Bonnivard et al. 2015). 

Light profile: We use a generalised Zhao-Hernquist profile 
(Hemquist 1990; Zhao 1996) for the stellar number density: 

Zhao / \ _ 

^ ” {r/r*)~t[l + {r/rt)°-]V^-~<)/°‘ ’ 

the five free parameters of which are the normalisation v*, the scale 
radius rJ, the inner slope 7 , the outer slope /3, and the transition 
slope a. Many studies have used less flexible parametrisations (e.g.. 
King, Plummer, or exponential profiles), but the use of these can 
bias the calculated astrophysical factors (Bonnivard et al. 2015). 


2.2 Likelihood functions 

2.2.1 Binned and unbinned analyses 


Before fitting the actual dSph kinematic data, we tested both a 
binned and an unbinned likelihood function on a set of mock 
data (mimicking ‘ultrafaint’ and ‘classical’ dSphs, see Appendix 
A). Both methods have been used in the literature, but to date, 
no systematic comparison has been undertaken to test the mer¬ 
its and limits of each approach (binned analyses can be found in 
Strigari et al. 2007; Charbonnier et al. 2011; unbinned in Strigari 
et al. 2008; Martinez et al. 2009; Geringer-Sameth, Koushiappas 
& Walker 2015). For the binned analysis, the velocity dispersion 
profiles (Jobs (7?) are built from the individual stellar velocities (see 
Section 3), and the likelihood function we use is: 


Wins ,_l/2 

/.bin -pr (27r) 

1=1 ^ ' 


1 ( <^obs(-^i) 

2 \ ISafiRi) 


, (8) 


where 

ai= aobs{Ri}v(^ [ap{Ri + ARi)—ap{Ri — ARi)]\. (9) 


The quantity Aaobs{Ri) is the error on the velocity dispersion at 
the radius Ri, and ARi is the standard deviation of the radii dis¬ 
tribution in the i-th bin. This likelihood allows the uncertainties on 
both (Jobs and R for each bin to be taken into account. 

For the unbinned analysis, we assume that the distribution of 
line-of-sight stellar velocities is Gaussian, centred on the mean stel¬ 
lar velocity v. The likelihood function reads (Strigari et al. 2008): 


£ 


unbin 


fJ^y/df{Ri)+^i 


exp 


U ivi-yf Y 
2\a^{R,) + AlJ 


, ( 10 ) 


where the dispersion of velocities at radius Ri of the i-th star comes 
from both the intrinsic dispersion ap{Ri) from equation (3) and the 
measurement uncertainty A„.. 

As detailed in Appendix A the unbinned analysis reduces the 
statistical uncertainties on the astrophysical factors, particularly for 
the ‘ultrafaint’ dSphs, without introducing biases. In the remain¬ 
der of the paper we therefore favour the unbinned analysis and the 
binned likelihood is used only to cross-check our results. 


2.2.2 Analysis with and without membership probabilities 


iBaes / \ _ do T fioo{r/Ca')^ 

'ani (f) z ' 7 7 r~ 

1 + (r/ra)-) 


Kinematic samples are often contaminated by interlopers from 
the Milky Way (MW) foreground stars. Different methods can be 
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used in order to remove those contaminants, based e.g. on sigma¬ 
clipping, virial theorem (Klimentowski et al. 2007; Wojtak & Lokas 
2007), or expectation maximisation (EM) algorithms (Walker et al. 
2009; Martinez et al. 2009). The latter allows in particular the es¬ 
timation of membership probabilities for each star of the object. 
These probabilities can be used as weights when building the veloc¬ 
ity dispersion profile (Tahs{R) in the binned case, or used directly 
in the unbinned likelihood, where equation (10) becomes 


i^unbin_ 

L-w - 



i = l 




exp 


if (vi — v)^ 
2[aliRi) + Al 


( 11 ) 


with Pi the membership probability of the i-th star. Another op¬ 
tion when dealing with foreground contamination is to use the un¬ 
weighted likelihoods analysis only with stars 

having large-enough membership probabilities (typically. Pi > 
0.95). Whenever membership probabilities are available, we test 
both methods (see Section 4.4). 


2.3 Astrophysical factor for annihilation and decay 

The 7 -ray differential flux from DM annihilation (resp. decay) in 
a dSph galaxy, measured within a solid angle Afl, is (Bergstrom, 
Ullio & Buckley 1998) 

^ X J(AD) (resp. x D(AT!)) . (12) 

The quantity (jfj(E-y) (resp. (j/j^(E^)'j is sensitive to the particle 
physics, e.g. the annihilation or decay channel. We focus here on 
the ‘astrophysical factor’, J (resp. D), 


J = JJ (resp. dldO^ , (13) 


which corresponds to the integration along the line-of-sight of 
the DM density squared (resp. DM density) over the solid angle 
AQ = 27 r X [1 — cos(aint)], Oint being the integration angle. 
A precise evaluation of this quantity is necessary for setting robust 
constraints on the properties of the DM particle, and is also a useful 
proxy to rank possible targets according to the magnitude of their 
flux. In order to compute the astrophysical factor, both the density 
distribution of the DM halo and its extent are required to be known. 
Results may be sensitive to the choice of the latter, as discussed 
further in Section 4.3 and Appendix B. 

N-body simulations within the context of ACDM cosmology 
have shown that DM halos should contain a large number of smaller 
sub-halos. Such sub-structures could significantly increase the J- 
factors, but the smaller the host halo mass, the less boosted is the 
signal (e.g., Sanchez-Conde & Prada 2014). For DM halos typical 
of dSphs, Charbonnier et al. (2011) found no significant impact of 
the sub-structures on the J-factors so we neglect their contribution 
here. 

J- and D-factors of dSphs were found to be best constrained 
at the so-called ‘critical’ integration angle (Walker et al. 2011; 
Charbonnier et al. 2011; Bonnivard et al. 2015). This is related to 
the half-light radius of the dSph, r^, and to its distance d, and dif¬ 
fers for J- and D-factors: ai ~ 2rhld, while ~ rh/d (see 
figure 4 for illustration). 

All calculations of astrophysical factors are done with the 
CLUMPY code (Charbonnier, Combet & Maurin 2012). A new 
module has been specifically developed to perform the Jeans anal¬ 
ysis, and this upgrade will be publicly available in the forthcoming 
second release of the software (Bonnivard et al., in preparation). 


The J- and D-factors obtained for the eight ‘classical’ and thirteen 
‘ultrafaint’ dSph galaxies are presented in Section 5. 


3 DWARF SPHEROIDAL GALAXY DATA 


3.1 Surface brightness data 


For each dSph, we estimate v{r) by fitting publicly available pho¬ 
tometric data with a Zhao-Hernquist model (equation 7), where r 
is the distance from the dwarf’s centre in 3D. Due to differences 
in the nature of the available data for ‘classical’ and ‘ultrafaint’ 
dSphs, we adopt different approaches for applying this model to 
each category. 

For ‘classical’ dSphs, the most homogeneous data sets for es¬ 
timating structural parameters remain those of Irwin & Hatzidim- 
itriou (1995, ‘IH95’ hereafter). IH95 tabulate stellar surface density 
profiles in terms of stars counted within concentric elliptical annuli, 
with each ellipse (with semi-major and semi-minor axes a and b 
respectively) having the same ellipticity and orientation that IH95 
estimate for the dSph as a whole. IH95 measure global ellipticities 
ranging from e = 1 — b/a — 0.13 (Leo II) to e = 0.56 (UMi), with 
a median of e = 0.32. Because our Jeans models assumes spherical 
symmetry, we transform IH95’s elliptical annuli into circular annuli 
by replacing the ‘elliptical radii’ in their data tables with geometric 
means, Rgm = '/ab. 

To the circularised, binned surface density profiles, S(J?), we 
then fit 2D projections of v{r) according to the likelihood function 


JVbi„= 

£i oc exp 

i=l 


1 {E{Ri) — Eniodel(7?i))^ 


(14) 


where is the Poisson error associated with the number of 

stars counted in the i**' bin and the model surface density is 

viv^T 

Smodel(7?) = 2 / =dr + Sbkd, (15) 

JR vf — R 

i.e., the sum of the projection of i/(r) and a uniform background 
density. 

For the ‘ultrafaint’ satellites, the largest homogeneous data 
set is from the Sloan Digital Sky Survey, which provides posi¬ 
tions, colours and magnitudes of individual stars detected as point 
sources. For each ‘ultrafaint’ satellite, we identify possible mem¬ 
bers as red giant branch (RGB) candidates, which we define as 
point sources whose g — r colours place them within 0.25 dex of 
the Dartmouth isochrone (Dotter et al. 2008) calculated for a stel¬ 
lar population with age 12 Gyr, metallicity corresponding to the 
mean value estimated from spectroscopy, and shifted by the dis¬ 
tance modulus estimated for that satellite (McConnachie 2012). To 
the unbinned distribution of projected positions for N RGB can¬ 
didates, we fit 2D projections of i>{r) according to the likelihood 
function 


N 

£2 OC T,model{Ri), (16) 

i=l 


where Eniodei(7?i) is defined as in equation (15). 

For both ‘classical’ and ‘ultrafaint’ satellites, we adopt uni¬ 
form priors on all model parameters and then use the software 
package MultiNest (Feroz & Hobson 2008; Feroz, Hobson & 
Bridges 2009; Feroz et al. 2013) to sample the posterior probability 
distribution (PDF) of parameter space. In a separate contribution. 
Walker (in preparation) will discuss detailed results from these fits. 
For our present purposes, we use these samples from the posterior 
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PDFs to propagate uncertainties in v{r) through our estimation of 
the dark matter density profiles (see Section 4). 


3.2 Kinematic data 

We use stellar-kinematic data, in the form of projected positions 
and line-of-sight velocities for individual stars, compiled from the 
literature. For all galaxies except Draco, we use the same data sets 
as analysed by Geringer-Sameth, Koushiappas & Walker (2015) 
(who provide a detailed description). For Draco, we adopt the data 
set of Walker, Olszewski & Mateo (2015), which includes measure¬ 
ments for ~ 500 members. 

The kinematic samples for the eight ‘classical’ dSphs include, 
for each star, a probability of membership. Pi, that is estimated 
using an expectation-maximisation algorithm (Walker et al. 2009). 
The sample for Segue I (Simon et al. 2011) includes membership 
probabilities that are estimated in two ways: via the same EM al¬ 
gorithm, and alternatively using a Bayesian analysis that considers 
the entire data set to be a mixture of Segue I and foreground pop¬ 
ulations. We use the latter estimation in our analysis - details of 
the Segue I case will be presented in a separate contribution (Bon- 
nivard, Maurin & Walker, in prep). Data for the remaining ‘ultra¬ 
faints’ do not include membership probabilities, but rather a bi¬ 
nary classification of each star as member or nonmember based on 
velocity and line-strength criteria (e.g., Simon & Geha 2007). We 
treat these classifications as membership probabilities in our analy¬ 
sis, but their values necessarily are either Pi (non-member) or 
Pi = 1 (member). 

In order to estimate velocity dispersion profiles, we divide 
each data set, consisting of Wmem = member stars, into 

~ V-^mem bins that each contain a number of stars whose mem¬ 
bership probabilities add to ~ 'JNme-m- For each bin, we estimate 
velocity dispersion using the maximum-likelihood procedure de¬ 
scribed by Walker et al. (2006), with membership probabilities in¬ 
troduced as weights on each star. Our results are based on the un¬ 
binned analysis (see Section 2.2), and we use these velocity disper¬ 
sion profiles only for cross-check purpose. 


4 ANALYSIS SETUP 
4.1 MCMC analysis 

For each dSph, we perform the Jeans analysis to find the DM den¬ 
sity profile parameters that best fit the stellar kinematics, and de¬ 
termine their uncertainties. The likelihood functions of the analysis 
(Section 2.2) have seven free parameters (see below). To efficiently 
explore this large parameter space, we use an MCMC technique, 
based on Bayesian parameter inference, which allows us to sample 
the PDF of a set of free parameters using Markov chains. To this 
purpose, we use the Grenoble Analysis Toolkit (GreAT) (Putze 
2011; Putze & Derome 2014), which relies on the Metropolis- 
Hastings algorithm (Metropolis et al. 1953; Hastings 1970). The 
posterior distributions are obtained after several post-processing 
steps (burn-in length removal, correlation length estimation, and 
thinning of the chains), allowing a selection of independent sam¬ 
ples, insensitive to the initial conditions. From these PDFs, credi¬ 
bility intervals (CIs) for any quantity of interest can be easily com¬ 
puted. 

For each step of the chains, we randomly select a light pro¬ 
file parametrisation from the accepted configurations of the MCMC 


Table 1. Range of uniform priors used for the DM density and velocity 
anisotropy profile parameters. Note that all models must satisfy the Global 
Density-Slope Anisotropy Inequality (equation 17), which reduces the ef¬ 
fective range of the anisotropy parameters. See text for details. 


Quantity Profile 

Parameter 

Prior range 

DM density ‘Einasto’ 

logio(P- 2 /M 0 kpc-3) 

[5,13] 

equation (5) 

logic (I'-z/kpc) 

[logic «).l] 


a 

[0.12,1] 

Anisotropy ‘Baes & van Hese’ 

Po 

[-9,1] 

equation ( 6 ) 

fioa 

[-9,1] 


logic 

[-3,1] 


7 

[0.1,4] 


analysis done previously for the surface brightness data (see Sec¬ 
tion 3), and use it for computing the velocity dispersion ap (equa¬ 
tion 3). This effectively propagates the light profile uncertainties to 
the posterior distributions of the DM and anisotropy parameters. 


4.2 Free parameters and priors 


The free parameters of the analysis are the three parameters of the 
Einasto DM density profile, and the four parameters of the Baes 
& van Hese (2007) velocity anisotropy profile, for which we adopt 
uniform priors. We follow Bonnivard et al. (2015) who after a thor¬ 
ough study of mock data identified an optimal prior combination 
that mitigates several biases introduced by the Jeans analysis. Ta¬ 
ble 1 summarises the ranges used on each parameter’s prior. 

DM density profile. Bonnivard et al. (2015) suggested two 
cuts on the DM density profile priors in order to tighten the con¬ 
straints on the astrophysical factors without introducing bias. Eirst, 
the scale radius r _2 is forced to be larger than the scale radius of 
the stellar component, r*. This drastically reduces the upper CIs 
on the astrophysical factors for ‘ultrafaint’ dSphs. A second cut on 
the logarithmic slope a, a 0.12, is also advocated in this previ¬ 
ous work. We have chosen to use both cuts in our analysis, so as to 
obtain the most stringent and robust estimates of the astrophysical 
factors. 

Velocity anisotropy profile. The priors we use for the 4 param¬ 
eters of the Baes & van Hese (2007) anisotropy profiles are also 
those of Bonnivard et al. (2015). The generality of this parametri¬ 
sation avoids the bias from the use of more specific anisotropy 
profiles, especially for large data samples. We also implement the 
Global Density-Slope Anisotropy Inequality (Ciotti & Morganti 
2010), which ensures that solutions to the Jeans equation given by 
our MCMC analysis correspond to physical models with positive 
phase-space distribution function: 


/3ani(r) ^ 


1 dlog v{r) 

2 dlog(r) 


(17) 


This reduces the range of allowed velocity anisotropy parameters 
depending on the stellar number density o{r). 


4.3 Size of the DM halo 

The extension of the DM halo is needed when computing the J- (or 
D-) factor (equation 13). The latter reaches a maximum when the 
integration angle Oint corresponds to the angular size of the halo 
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Figure 1. Velocity dispersion profiles cjp of the eight ‘classical’ dSphs: data (symbols) and reconstructed median and 95% Cls (blue lines). These profiles are 
shown for illustration; our results are based on an unbinned analysis. 


and saturates beyond (see figure 4, where the J-factor obtained for 
Fornax is plotted as a function of the integration angle Oint; the 
median value is seen to saturate above ~ 1°). 

There is no clear criterion to define the size of DM halos 
hence we adopt two different approaches for each dSph galaxy and 
for each set of DM parameters accepted by the MCMC. The first 
method considers the tidal radius n to be a good estimator of the 
halo size (as shown by N-body or hydrodynamical simulations — 
see e.g. Springel et al. 2008; Mollitor, Nezri & Teyssier 2015); this 
is computed as: 


Mhalo(n) 


[2 — dinMiviw/dInr(d)] x MMw(d) 


(1/3) 


xd, (18) 


where Mmw (d) is the mass of the MW enclosed within the galac- 
tocentric distance d of the dSph and Mhaio is the mass of the dSph 
galaxy. A second method to estimate the size of a DM halo consists 
in determining the radius req where the halo density is equal to the 
density of the MW halo, namely, 

PDM(req)=P^,S'(d-req). (19) 

We have used both a Navarro-Frenk-White (Navarro. Frenk & 
White 1997; Battaglia et al. 2005) and an Einasto profile (Navan'o 
et al. 2004; Springel et al. 2008) for the MW density, with no im¬ 
pact on the results. We also find the above two estimates of the dSph 
galaxy size to be comparable, leading to very similar astrophysi- 
cal factors (see Appendix B). Note that other assumptions can also 
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Figure 2. Velocity dispersion profiles o-p of the thirteen ‘ultrafaint’ dSphs: 
are shown for illustration; our results are based on an unbinned analysis. 

be made to estimate the dSph halo size. For instance, in order to 
be conservative, Geringer-Sameth, Koushiappas & Walker (2015) 
used the outermost observed star as truncation radius for comput¬ 
ing the astrophysical factors. However, this can underestimate the 
J- and D-factors, and may moreover underestimate the credibility 
intervals (see Appendix B). 

4.4 Membership probabilities and impact of contamination 

J- and D-factor reconstruction through Jeans analysis has been 
previously studied and optimised for contamination-free mock data 
(Charbonnier et al. 2011; Bonnivard et al. 2015). However, actual 
observations yield kinematic samples that may be contaminated by 


(symbols) and reconstructed median and 95% CIs (blue lines). These profiles 


field stars belonging to the Milky Way or to a Galactic stream. The 
conventional approach to handle these interlopers relies on the ar¬ 
bitrary definition of some threshold separating members from out¬ 
liers (sigma-clipping method, see e.g. Yahil & Vidal 1977). Expec¬ 
tation maximisation (EM) algorithms (Walker et al. 2009) differ 
from sigma-clipping methods as they provide membership proba¬ 
bilities Pi for each star of the sample, which can be used as weights 
in subsequent analyses, e.g., equation (11). The EM algorithm was 
shown to provide accurate and reliable membership probabilities in 
most cases, although some failures may occur for samples present¬ 
ing the heaviest contamination and the most overlapping velocity 
distributions (Walker et al. 2009). 

In order to investigate whether residual contamination affects 
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Figure 3. J- and _D-factors at the critical integration angle (top and bottom 
panel respectively) for dSph galaxies whose Pi values are available. Shown 
are values calculated using equation (11), i.e. with membership probabilities 
as weights in the Jeans analysis (red triangles), or using all stars passing the 
cut Pi > 0.95 (blue squares). 


the J-factor values, we either use the star membership probabil¬ 
ity as weights in the likelihood function equation ( 11 ), or only re¬ 
tain stars which almost certainly belong to the dSph (Pi > 0.95). 
Figure 3 compares these two approaches for dSphs with Pi val¬ 
ues available from the literature (eight ‘classical’ and Segue I). The 
same J-factors (D-factors) are reconstructed in both cases, except 
for Fornax and Segue L For these two objects, the discrepancy hints 
at the presence of contamination in their samples that is not cap¬ 
tured by the Pi indicator. The case of Segue I will be thoroughly 
discussed in Bonnivard, Maurin & Walker (in prep.), and we refer 
the reader to Appendix C for the ‘classical’ dSphs. Investigation of 
the stellar contamination issue is done by comparing reconstructed 
J-factors to their tme values, for a set of mock data presenting dif¬ 
ferent levels of Milky Way or stream contamination. Doing so, we 
find that: 

• J-factors can be robustly reconstructed when large enough 
samples of stars are available (‘classical’ dSphs), whereas small 
samples (‘ultrafaint’ dSphs) are more sensitive to contamination, 
with their J- and JJ-factors more likely to overshoot the true value; 

• discrepant J- or D-factors from the P^-weighted and 

cut analyses hint at high levels of contamination. The PJ° ®®-cut 
analysis is found to give generally more conservative results (large 
CIs, but encompassing the true value) while the Pi-weighted anal¬ 
ysis tends to overshoot (small CIs, true value outside CIs). 

Therefore, in the remainder of the paper our results are based 
only on stars with Pi > 0.95, whenever this information is avail¬ 
able. A direct consequence of this is that Segue I (among the most 


favoured targets) may become one of the least reliable targets to set 
constraints on DM (Bonnivard, Maurin & Walker in prep.). Unsur¬ 
prisingly this confirms that the ‘classical’ dSphs are the most robust 
targets as they are little affected by contamination. 


5 RESULTS 

Using the MCMC analysis described in the previous Section, we 
fit the velocity data of the 8 ‘classical’ and of 13 ‘ultrafaint’ dSph 
galaxies with the 7 parameters (3 parameters for the Einasto DM 
profile and 4 for the Baes & Van Hese velocity anisotropy) required 
in our Jeans analysis. For each point in the chains, any relevant 
quantity may be computed and its median value and credibility in¬ 
tervals estimated from the resulting distribution. 

This is true of the reconstructed velocity dispersion profile of 
each dSph, the median and 95% CIs of which are plotted in solid 
and dashed blue lines in figures 1 and 2d The reconstructed profiles 
and CIs appear to always provide a good representation of the data, 
with much wider CIs for ‘ultrafaint’ dSphs (up to a factor ~ 5 at 
large radii, compared to a factor ~ 2 for ‘classical’ dSphs) because 
of the sparsity of their stellar data. 

The three Einasto profile parameters are used to compute the 
DM annihilation J-factors and decaying DM D-factors using equa¬ 
tion (13). Figure 4 displays the median value and 95% CIs of the 
J- and U-factors, for the ‘classical’ dSph Fornax, as a function of 
Qint. Once the integration angle encompasses the whole halo, the 
astrophysical factors saturate. This figure shows the optimal inte¬ 
gration angle for which the CIs for the J- and D-factors are the 
smallest, i.e., where these astrophysical factors may be robustly de¬ 
termined despite our inability to constrain the inner slope of their 
DM profile. For Fornax, ~ 0.28° as may be seen from the 
figure. ASCII files of the median values, 68 % and 95% CIs of 
J(Q:int) and D(aint) for all the dSphs discussed in this paper may 
be retrieved from the Supporting Information submitted with this 
paper. 


5.1 J and D-factors of dSphs vs Galactic background 

All the dSph satellite galaxies of the MW are embedded in its DM 
halo. Hence in both the annihilating and decaying DM scenario, our 
Galaxy’s DM halo will provide a background signal of the same na¬ 
ture as that of the targeted dSph galaxy. This consideration is quite 
independent of any diffuse 7 -ray emission of astrophysical origin 
which we do not discuss here. For simplicity, we also ignore the 
extragalactic signal originating from DM annihilations or decays 
on cosmological scales (and integrated over all redshifts). 

The MW halo’s astrophysical J- and D-factors are computed 
with the Clumpy code, assuming the following characteristics: 

• We use an Einasto profile to model the smooth DM distribu¬ 
tion, scaled to the local DM density (p© = 0.3 GeV cm“®); 

• We include the contribution of a population of DM clumps, 
having a cored spatial distribution and a mass distribution 
dN/dM oc M~^'^ (Springel et al. 2008) (amounting to ~ 10% 
of the MW’s total mass); 


^ For display purposes, the binned analysis has been used in these figures, 
i.e. binned data and likelihood function. We have checked for each 
dSph that the results obtained using the binned or unbinned analysis are 
consistant. 
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Figure 4. J- and D-factor (top and bottom panel respectively) as a function 
of the integration angle aim for the ‘classical’ dSph Fornax. 

• We assume that the mass-concentration relation® has a log¬ 
normal distrihution (Sanchez-Conde & Prada 2014). 

Figure 5 displays the J- (left) and D-factors (right) of the 
21 dSph galaxies (symbols) studied in this paper as a function 
of their angular distance from the Galactic centre. The solid blue 
line corresponds to our estimate of the contribution from the MW 
DM halo. This has been repeated for three integration angles, 
Q!int = 0.01°,0.1°, 1°. This figure clearly illustrates the loss of 
contrast between the dSph target and the MW background as the 
integration angle is increased. Indeed, the background (exotic or 
not) is oc Omti this is not so for the dSphs where the astrophysical 
factor is very centrally peaked, especially for J. For the J-factor, 
most of the dSphs appear an order of magnitude or more above the 
background for aint = 0 . 01 °, while this is hue only of a couple of 
them for aint = 1°. For decay (right column), even for small inte¬ 
gration angles, the contrast is always smaller than that of annihila¬ 
tion. Therefore, for large integration angles which may be dictated 
by instrumental resolution, it would be a better strategy to look di¬ 
rectly for the MW’s DM halo signature, rather than at a specific 
dSph galaxy. 


® The mass-concentration relation is a fundamental ingredient to relate a 
halo mass and size to the scale density and radius required by the DM den¬ 
sity parametrisation. 


5.2 J-factor: ranking of the dSphs and comparison to other 
works 

Putting aside the notion of contrast mentioned above, the astrophys¬ 
ical factors are the relevant proxies to determine whether or not a 
given dSph is potentially interesting for indirect detection. Table 2 
summarises our findings and we compare these results in figure 6 
to other studies in the literature. 

Ranking at aint = ~ 2rhld. The top panel in figure 6 shows 

the J-factors ordered according to their median from the faintest 
(Leo5) to the brightest (UMa2), when integrating the signal up to 
the optimal angle aint = ctc, where the J-factors have the smallest 
error bars. To account for possible systematics from triaxiality of 
the DM halo (which depends on the l.o.s. orientation of each dSph, 
see Bonnivard et al. 2015), the error bars (’[]’ symbols) combine 
a 0.4 dex uncertainty in quadrature with the 68 % CIs. Our up¬ 
dated analysis (compared to the results presented in Charbonnier 
et al. 2011) still prefers UMi among all the ‘classical’ dSph while 
Coma and Ursa Major 2 are the most promising ‘ultrafaint’ dSphs 
(Segue I falls among the least interesting targets). For each dSph 
galaxy, the optimal angle Q;)f is quoted above the data point, while 
the distance to the dSph is quoted below. Unsurprisingly, the most 
promising galaxies are also among the closest. The estimated un¬ 
certainties on this plot provide a clear signature of our data-driven 
approach: the ‘ultrafaint’ dSphs have significantly larger error bars 
than their ‘classical’ counterparts. 

Ranking at Qint = 0.5° and comparison to other works. The 

bottom panel in figure 6 compares our results to other existing stud¬ 
ies for a fixed integration angle aint = 0.5°, typical of the Fermi- 
LAT angular resolution in the GeV range. First, comparing the top 
and bottom panels shows no drastic ordering changes among the 
dSph galaxies but for a few inversions. Second, the overall trend ap¬ 
pears to be preserved between the different studies, with the same 
objects having the highest J values. Nonetheless, a closer inspec¬ 
tion shows significant differences both in the values of the J-factors 
and the size of the error bars. The differences observed between our 
analysis and others is understood as follows: 

• In Charbonnier et al. (2011), we conducted the same exer¬ 
cise on the 8 ‘classical’ dSph galaxies (green circles), but with a 
more constrained Jeans analysis (using a constant anisotropy pro¬ 
file and a Plummer light profile). UMi was found to be the most 
promising target among the ‘classical’ dSph, while Leo 2 had the 
highest median J{ad) but with much larger uncertainties. For all 
the objects but Leo 2 (which is very uncertain anyway), the Char¬ 
bonnier et al. ( 2011 ) values are slightly lower than their updated 
version. Changing the anisotropy and light profiles to more flexi¬ 
ble parametrisations in the current analysis is the main reason for 
the differences between this and our previous study. Overall, we 
find larger J-factors whenever the light data require an outer slope 
steeper than that given by the Plummer profile used in Charbonnier 
et al. (2011). Note that this effect was already observed on mock 
data in Bonnivard et al. (2015). The most striking example is the 
case of Draco. 

• The Fermi collaboration (Fermi-LAT Collaboration 2014, 
2015) reported the 7 -ray observations of 25 dSph galaxies and 
conducted a stacking analysis of 15 of them to set constraints on 
the DM annihilation cross section. These authors provide J-factors 
of 18 dSph galaxies overlapping with our sample (orange trian¬ 
gles), obtained using a two-level Bayesian hierarchical modelling 
that constraints the entire population of MW dSphs simultaneously 











10 Bonnivard, Combet, Daniel et al. 



Figure 5. Annihilation (left) and decay (right) factors of the 21 dSphs studied (symbols) versus the Galactic DM halo background (blue line), as function of 
the angular distance from the Galactic centre. Top to bottom panels correspond to three integration angles: 0.01°, 0.1°, and 1°. See Section 5.1 for details. 


(Martinez 2013). Their values do not show any particular trend 
when compared to ours, with most CIs overlapping. The most strik¬ 
ing difference concerns the size of their error bars, which remain 
roughly the same regardless of the nature of the object (‘classical’ 
or ‘ultrafaint’). This is very likely related to their two-level hier¬ 
archical analysis, where the entire population of MW dSph galax¬ 
ies is to some extent assumed to share the same properties. The 
constraints coming from the ‘classical’ dSphs are therefore redis¬ 
tributed to the ‘ultrafaint’ dSphs. Note also that a NFW profile for 
the DM density was assumed, but the analysis found to be fairly 
insensitive to this choice. 

• In Geringer-Sameth, Koushiappas & Walker (2015), the au¬ 
thors provide J- and D-factors of 20 dSph galaxies obtained us¬ 
ing a data-driven Jeans analysis quite similar to ours. In particu¬ 
lar, they use the same unbinned likelihood as in equation (10), but 
with a Zhao DM density profile, a constant anisotropy and a Plum¬ 
mer light profile. However, they select the radius of the outermost 
observed star as truncation radius for computing the astrophysi- 
cal factors, which results in lower values (red triangles) than ours 
(blue squares). Segue II (Seg2), Hercules (Her), and Ursa Major I 
(UMal) are particularly affected by this choice, with 7?max/rh val¬ 
ues between ~ 2 and ~ 4, while we obtain J?max/?'h values > 10 
with our estimations of the halo size (see Section 4.3). Note that the 


underestimation of the halo size can also lead to an underestimation 
of the CIs (Appendix B). This partially explains the larger error 
bars we find, which are also a consequence of the more flexible 
anisotropy and light profiles parametrisations we use (Bonnivard 
etal. 2015). 


In summary, figure 6 (bottom) highlights the fact that the val¬ 
ues found for J-factors and their CIs can depend strongly on the 
underlying assumptions. We believe our present work to better re¬ 
flect the various sources of uncertainties that affect J-factor estima¬ 
tions, following the thorough validation of the method initiated in 
Bonnivard et al. (2015) and concluded in Appendices A and B. Our 
data-driven analysis naturally implies larger error bars for dSphs 
with fewer stars, and vice-versa. Also, thanks to a more flexible 
parametrisation of light profile, we find higher J-factors for some 
targets, which may lower further the upper limits for the DM an¬ 
nihilation cross section (Fermi-LAT Collaboration 2014; Geringer- 
Sameth, Koushiappas & Walker 2014; Fermi-LAT Collaboration 
2015). 
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from triaxiality of the dSph galaxies (Bonnivard et al. 2015). Bottom', comparison of the J-factors to other works, with aint = 0.5°. See also Section 5.4 for a 
critical discussion of the targets most favoured by our analysis. 


5.3 D-factor: ranking of the dSphs and comparison to other 
works 

Dark matter decay is less often considered than annihilation, how¬ 
ever recent observations of an unidentified X-ray line at 3.55 keV in 
galaxy clusters has generated increasing interest in this possibility 
(e.g., Bulbul et al. 2014a; Boyarsky et al. 2014b). 

Ranking. The blue squares in figure 7 and the three rightmost 
columns of table 2 give an overview of the D-factors computed 
here. First, comparing the top panels of figures 6 and 7, we find that 
the ordering of the most promising targets changes significantly 
whether focusing on DM annihilation or decay, even though Ursa 
Major II remains the best candidate for Oint = Q?. Furthermore, 
the two panels in figure 7 show that changing the integration angle 
for a decaying DM signal also has a strong impact on the ranking 


and on the error bars, more strongly than in the case of DM anni¬ 
hilation. In particular, for aint =0.1° (bottom panel), most targets 
have very similar D-factors and the increased error bars make the 
ranking less obvious. 

Comparison to other works. The availability of independently- 
derived D-factors for dSphs in the literature remains limited, mak¬ 
ing comparison less straightforward than in the case of annihilation. 

• Although not published in the Charbonnier et al. (2011) study 
which focused on J-factors only, the D-factors for the eight ‘clas¬ 
sical’ dSphs were also obtained from our original analysis setup. 
As in the case of annihilation, these values (green dots in figure 7) 
are systematically lower than that obtained by the present analysis 
and this is connected, as for J, to the choice of the light profile. 

• We also compare our results to those of Geringer-Sameth, 
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Table 2. Summary of results for the 21 dSph galaxies presented in this paper. The dSphs are ordered by distance and the columns correspond to (from left 
to right): name, distance, optimal angle for annihilation (a^ Ri ^ 2r(,/d), median J-factors and 68% (95%) CIs for aint = 0.01, 0.5, median 
D-factors and 68% (95%) CIs for = 0.01, 0.5, a^?. Note that a systematic uncertainty of ±0.4 (resp. ±0.3) must be allowed in order to the reflect the 
possible triviality of the dSph galaxies (Bonnivard et al. 2015). 
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„ +0.5(+0.9) 
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Ursa Minor 

66 

0.49 
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11 q+0.6(+1.1) 
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79 

0.38 

1 n 
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-1^ +0.1(+0.3) 
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86 

0.91 
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0.27 
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+0.2(+0.5) 

-*--‘--'J_o.i(_o.2) 
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Leo IV 

160 

0.08 
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0.08 
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0.11 
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0.05 
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.. .+0.4(+0.8) 
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Koushiappas & Walker (2015), noting again that the conservative 
choice made for the size of the DM halo in that study leads to a 
deficit in the D-factors compared to our values. This deficit is more 
apparent here than for annihilation as the outer regions of the DM 
halo contribute more to the D-factors than to the J-factors. 

• Following the claims of X-ray line detections in galaxy clus¬ 
ters, Malyshev, Neronov & Eckert (2014) looked for such a sig¬ 
nal in the dSph galaxies available in XMM-Newton data. In the 
absence of a signal, these authors used the mass derived from 
Wolf et al. (2010) ® and Geringer-Sameth, Koushiappas & Walker 


® Wolf et al. (2010) solve the spherical Jeans equation coupled to a MCMC 
technique, using a similar approach to ours but different profile (mass, light, 
anisotropy) parametrisations, to provide a robust mass estimate of several 
dSph galaxies. 


(2015) to set constraints on a sterile neutrino DM scenario. In do¬ 
ing so, rather than writing the D-factor as given by equation 13, 
these authors use the point-like approximation and use instead 
where d is the distance to the galaxy 
and is the enclosed mass.^ Wolf et al. (2010) provide the 

mass Mi /2 (and the corresponding error bars) contained within the 
deprojected half-light radius ri/ 2 , corresponding to an integration 
angle Qi /2 = t'i/ 2 /d- This angle is closely related to our defini¬ 
tion of the optimal integration for decay = rh/d, where rh ~ 


^ Malyshev, Neronov & Eckert (2014) do not explicitly call this quantity 
the D-factor; it simply appears as pail of the overall flux definition. In the 
point like approximation, this mass is the mass enclosed in a sphere of ra¬ 
dius «int X d. It is not the mass contained in the volume defined by the 
intersection of the line-of-sight cone and the dSph spherical halo. 
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Figure 7. Top\ D-factors and 68% CIs for aim = oiP and comparisons to other works Bottom: same but for oim = 0.1°: the ’[]’ symbols combine the 
68% statistical uncertainties and possible systematics (±0.3) from triaxiality of the dSph galaxies (Bonnivard et al. 2015). See also Section 5.4 for a critical 
discussion of the targets most favoured by our analysis. 


0.75 X r-i /2 is the projected half-light radius.® For each dSph in the 
Wolf et al. (2010) sample, we compute = Afi/ 2 /rf^ 

(orange triangles in figure 7). Despite the fact that ai /2 > 
our D-factors are generally higher than the ones derived from the 
Wolf et al. (2010) data. Computing M 1/2 from our MCMC chains 
provides, in general, values that are compatible with those of Wolf 
et al. (2010); this excludes the mass reconstruction from being the 
sole origin of the differences. The main remaining difference may 
lie in the point-like approximation, which does not account prop¬ 
erly for the full volume of the DM halo being intercepted by the 
line-of-sight. While this may not be an inappropriate assumption 


® Wolf et al. (2010) provides useful fitting formulae to relate the projected 
and deprojected half-light radii, for a variety of light profiles. 


for the strongly peaked annihilation signal, it results in a significant 
deficit for the D-factors. 


5.4 Discussion 

In the previous sub-sections, we presented a ranking of the Milky 
Way dSph satellites as potential targets for dark matter annihila¬ 
tion/decay surveys based on the estimated values and uncertainties 
of their J-factors or D-factors. However, it is important to note that 
while our analysis has marginalised over many of the modelling un¬ 
certainties, issues such as the dynamical status of individual dSphs 
and evidence pointing to cored profiles in a number of the ‘classi¬ 
cal’ dSphs are not accounted for. 

First, it is possible that some of the dSphs, especially the ‘ul¬ 
trafaints’, are not currently in dynamical equilibrium. In particular. 
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a number of authors have presented evidence that the UMa2 dSph, 
which occupies the top position in both the J- and D-factor rank¬ 
ings, is currently experiencing strong tidal disturbance by the Milky 
Way (Fellhauer et al. 2007; Munoz, Geha & Willman 2010; Smith 
et al. 2013). While it is not clear that this has inflated the veloc¬ 
ity dispersion of UMa2, some caution is advisable before selecting 
UMa2 as a prime candidate for indirect dark matter detection sur¬ 
veys. 

Secondly, the precise nature of some ‘ultrafaint’ galaxies is 
still uncertain—it is possible that some are more closely related to 
star clusters and do not, in fact, contain dark matter. For example, 
the most recent study of Will describes it as ”A Probable Dwarf 
Galaxy with an Irregular Kinematic Distribution” and notes that 
foreground contamination and unusual kinematics make the deter¬ 
mination of its dark matter content difficult (Willman et al. 2011). 
Target selection for dark matter surveys must take account of such 
uncertainty when weighting candidates for study. 

The Sextans, Ursa Minor and Draco ‘classical’ dSphs feature 
in the top five of both the J- and D-factor rankings. The larger kine¬ 
matic samples in these objects make the Jeans modelling more ro¬ 
bust (as indicated by their smaller confidence intervals in figures 6 
and 7). However, in the case of Sextans (Lora et al. 2013; Kleyna 
et al. 2004) and Ursa Minor (Pace et al. 2014; Kleyna et al. 2005), 
it has been suggested that the presence of kinematic suhstructures 
indicate that the dark matter halos of these objects are cored rather 
than cusped. In principle, this could be taken into account via priors 
on the halo slope within our Bayesian analysis, but as shown in fig¬ 
ure 15 of Charbonnier et al. (2011 ), this would not actually change 
the conclusions for most integration angles (i.e., aint > 0.1°). 

Thus, of the seven dSphs which are within the top five of the J 
and D rankings, additional evidence for 4 of them suggests that the 
results of equilibrium dynamical modelling might not be sufficient 
to characterise their suitability as targets for indirect dark matter 
detection. Clearly, it is imperative that further data are obtained on 
potential candidates before future surveys select their target lists. 
New non-equilibrium modelling approaches such as that presented 
in Ural et al. (2015, in press) will also be important in placing the 
target selection on a secure footing. 


6 CONCLUSIONS 

Dwarf spheroidal galaxies have been widely targetted in for 
searches for annihilating dark matter in the Galaxy. This has en¬ 
abled 7 -ray telescopes to set very stringent limits on the DM anni¬ 
hilation cross section which are now beginning to impinge on the 
cosmologically preferred value for a thermal relic. Reliable esti¬ 
mates of the dSph J-factors and associated error budgets are clearly 
crucial in this regard. This is especially true for ‘stacking’ analyses 
that use data on several dSph galaxies simultaneously to improve 
the sensitivity. Ranking of the dSphs, according to their J- (and D-) 
factors, is also mandatory to optimise the strategy of pointed obser¬ 
vations. In case of a positive detection in a given target, this will 
inform the strategy for subsequent observations aiming to validate 
the DM hypothesis. 

This study — follow-up of our previous effort (Charbonnier 
et al. 2011) — extends and improves the reconstruction of the as- 
trophysical factor for dSph galaxies in several ways: 

• We use the optimised analysis setup proposed in Bonnivard 
et al. (2015): the parametrisation of the ingredients of the Jeans 
analysis are kept as general as possible to minimise biases. In the 


spirit of Charbonnier et al. (2011), we adopt very weak priors to 
have as data-driven an analysis as possible. 

• We rely on an improved analysis of light profiles (Walker et 
al., in prep.) for better J-factor reconstruction (Bonnivard et al. 
2015), and also include recent kinematic data (e.g. for Draco, see 
Walker, Olszewski & Mateo 2015) in the analysis. 

• We test the impact of the choice of the likelihood function 
on the results: the performance and consistency of binned and un¬ 
binned analyses are validated on real and mock data. Furthermore, 
contamination from foreground stars and any associated impact on 
the J-factors are also investigated, using the membership probabil¬ 
ity of stars when available. 

• In addition to the 8 ‘classical’ dSph galaxies in Charbonnier 
et al. (2011), re-analysed here, the results for 13 ‘ultrafaint’ dSphs 
are now also provided — both the J-factors for annihilating DM 
and the D-factors for decaying candidates. 

The most important result of our study is the ranking (median 
and CIs) of the astrophysical factors shown in figures 6 (annihila¬ 
tion) and 7 (decay), and summarised in Table 2. ASCII files for 
the median and 68% and 95% CIs for a large range of integration 
angles can be retrieved from the Supporting Information submitted 
with this paper. Our findings can be summarised thus: 

(i) The unbinned Jeans analysis of stars whose membership 
probability is > 0.95 gives the most stringent constraints and is 
appropriate to deal conservatively with possible contaminations. 

(ii) Our J- and D-factors are in general consistent with other 
calculations, though several differences are observed: using a more 
flexible light profile (cf. the usually adopted Plummer profile) 
slightly increases the astrophysical factor for several dSphs; using 
a more flexible anisotropy profile (w.r.t. the usually assumed con¬ 
stant) slightly enlarges the J-factor CIs, providing more realistic 
uncertainties. This also mitigates possible biases (Bonnivard et al. 
2015). 

(iii) Uncertainties on the astrophysical factors (J and D) from 
this data-driven analysis are directly related to the sample size used 
for the analysis (large error bars for ‘ultrafaints’, small error bars 
for ‘classicals’). We believe this better accommodates a possible 
non-universality in the properties of these objects. 

(iv) The ranking of the targets (according to their median val¬ 
ues) slightly depends on the integration angle. At the optimal angle 
ai ~ 2rh / d, the ‘classical’ dSphs UMi and Draco are confirmed 
as the potentially-brightest and most favoured targets in terms of J- 
factors. The ‘ultrafaint’ objects UMa 2 and Coma outrank them, but 
suffer from larger uncertainties, and in particular their lower 95% 
CIs are lower than those of UMi or Draco. For decaying DM, Sex¬ 
tans appears as the brightest ‘classical’ dSph at af « r^/d, while 
UMa 2 remains the brightest ‘ultrafaint’ target. Not discussed here 
is the frequency-dependent astrophysical background that may af¬ 
fect this ranking. This would require an instrument-specific analy¬ 
sis that goes beyond the scope of this paper. 

(v) The astrophysical factors of Segue I might be highly uncer¬ 
tain due to probable stellar contamination and few kinematic data, 
and this object does not make it among the top ten targets (Bon¬ 
nivard, Maurin & Walker, in prep.). (Re)-analyses of membership 
probabilities for other ‘ultrafaints’ would be helpful to probe the 
level of contamination in these objects. 

The 9 new potential dSphs discovered in the DBS survey 
(DBS Collaboration et al. 2015; Koposov et al. 2015) and already 
searched for in Fermi-LAT data (Fermi-LAT Collaboration et al. 
2015; Geringer-Sameth et al. 2015) call for a continued effort on 
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this topic (see also Laevens et al. 2015; Martin et al. 2015; Kim 
et al. 2015 for three other recently discovered dSphs). Careful stud¬ 
ies of their astrophysical factors will likely be difficult due to the 
small amount of kinematic data expected for these objects. As un¬ 
derlined in this study, it is all the more important to explore the 
limitations of kinematic analyses for such small stellar samples. 
This may prove crucial to understand whether the recently pub¬ 
lished limits from 6 years of Fermi-LAT observations on 15 dSphs 
(Geringer-Sameth, Koushiappas & Walker 2014; Fermi-LAT Col¬ 
laboration 2015) can be significantly improved or not, and/or if bet¬ 
ter targets exist for the forthcoming Cherenkov Telescope Array 
(Actis et al. 2011). 
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APPENDIX A: UNBINNED ANALYSIS - VALIDATION ON 
MOCK DATA 

Bonnivard et al. (2015) thoroughly tested the Jeans analysis using 
binned velocity dispersion profiles (see equation 8), but unbinned 
analyses (equation 10) are also often used in the literature (e.g. Stri- 
gari et al. 2008; Martinez et al. 2009; Geringer-Sameth, Koushiap¬ 
pas & Walker 2015). Here, we conduct a careful comparison using 
a large set of mock data to determine the merits and limits of each 
approach, and select the optimal setup for our analysis. 

Mock data set To this purpose, we employed the same suite of 
mock data set used previously by Walker et al. (2011); Charbonnier 
et al. (2011); Bonnivard et al. (2015). It consists of 64 models cov¬ 
ering a large variety of DM density profiles (from cored to NFW- 
like cuspy profiles), stellar light profiles and velocity anisotropy 
values (between /Jani = —0.45 and /3ani = -1-0.3, with constant 
anisotropy). We refer the reader to the papers quoted above for a 
more complete description of this mock data set. For each model, 
we draw samples of W = 30 (small), 1000 (medium) and 10000 
(large) stars, mimicking ‘ultrafaint’, ‘classical’ and ‘ideally ob¬ 
served’ dSphs. No foreground contamination is added to the data 
sets, and all the objects are fixed to a distance d = 100 kpc. 

Analysis setup For each 64 mock models and sample sizes com¬ 
bination, we run the Jeans analysis using either the binned® or un¬ 
binned likelihood functions described in Section 2.2 (see Eqs. 8 
and 10 respectively). We fix the light and anisotropy parameters to 
their true values (i.e. the ones used to generate the mock data), in 
order to disentangle the effects of the likelihood functions to that 
originating from the physical parameters. Eor each model we then 
compute the J- and D-factors as a function of the integration angle 

® Note that we did not take into account the radius uncertainty on the ve¬ 
locity dispersion profiles {ARi in equation 9) in the Bonnivard et al. (2015) 
analysis. 



Figure Al. Median value (solid) and 95% CIs (dashed) of the J-factors as 
a function of the integration angle aint, reconstructed using either a binned 
(blue tilled circles) or an unbinned (red empty circles) Jeans analysis on a 
mock ‘ultrafaint’ dSph. Both analyses lead to very similar results, with the 
true J-factor (black) being encompassed within the reconstructed CIs. 

Effects of unbinned analysis Eirst, we find that using either the 
binned or the unbinned likelihood function leads to very similar 
results for the 64 models, regardless of the sample size. This is 
illustrated in figure Al, by comparing the J-factors obtained using 
either the binned (blue filled circles) or the unbinned (red empty 
circles) analysis to the true value (solid back line, computed using 
the DM parameters used to generate the mock data), for a typical 
‘ultrafaint’-like mock dSph. 

The main difference between the two analyses lies in the CIs, 
as the unbinned analysis is found to be more constraining than the 
binned approach for ‘ultrafaint-like’ dSphs. We show in figure A2 
the distributions of jmed^ ^ ^ lower 95% 

Cl to the median J-factor, obtained for the 64 mock ‘ultrafaint’ 
dSphs using either the binned (dashed blue) or the unbinned (solid 
red) analysis. The integration angle is set to 0.1 x Qc (top), etc 
(middle) and 10 x (bottom). Eor each integration angle, the 
mean value of j-95%ciy jmed educed by a factor two for the 
mock ‘ultrafaint’ dSphs when using the unbinned analysis. The ef¬ 
fect is much less pronounced for the medium (‘classical’-like) and 
large data samples. We note a similar effect on the upper 95 % CIs 
(not shown). Eor the D-factors, the effect is also present for the 
‘ultrafaint’-like dSphs, but less so, with a reduction of the CIs by 
~ 30% when using the unbinned analysis. 

Finally, we have checked that no bias is introduced when us¬ 
ing the unbinned analysis on the 64 mock models. We therefore 
advocate the use of the unbinned analysis when dealing with small 
data samples, as it allows a significant reduction of the statistical 
uncertainties. 


APPENDIX B: SIZE OE THE DM HALO 

As pointed out in Section 4.3, the way one defines the DM halo size 
could influence the values of the astrophysical factors. In order to 
quantify this effect, we performed several tests on the mock dataset 
described in Appendix A. 

Impact of halo size on J- and D-factors: tests on mock data 

First, using two typical mock models (with either cuspy or core 
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Figure A2. Distributions of j-95%Cly jmed^ obtained for the 64 mock 
‘ultrafaint’ dSphs using either the binned (dashed blue) or the unbinned 
(solid red) analysis. For the three integration angles considered («int = 
0.1 X Q!c, o-c and 10 X a^, from top to bottom), the unbinned analysis 
allows to reduce the mean J— 95%Cl^jmed ^ factor two. The effect is 
much less pronounced for larger data sets. 


DM density profile), and for three different choices of halo maxi¬ 
mum radius i?max, we compute the J-factor as a function of the in¬ 
tegration angle aint. Figure B1 compares the results to a reference 
value obtained with = 50 x r*, chosen arbitrarily. Unsur¬ 
prisingly it shows that an underestimation of the halo size leads to 
an underestimation of the J-factor, the effect being stronger for the 
core than for the cuspy DM profile. The J-factor can be underesti¬ 
mated by up to 70% at the critical angle Qc for a core profile and 
a halo size strongly underestimated (factor 25 too small). The ef¬ 
fect is even more important for the D-factors (underestimation by 
~ 80% at Qfc for the same model mentioned before, not shown). 

To evaluate the response of the CIs of the J-factors to the halo 
size, we run the Jeans analysis on the entire set of mock dSphs (see 
Appendix A). For each mock dSph, we fix three values of J7max, 
and compute the J-factors and their CIs from the reconstructed 
DM density profiles. We show in figure B2 the distributions of 
j+95%ci^jined three halo sizes, i.e. the ratio of the upper 

95% Cl to the median value, obtained for the 64 mock ‘ultrafaint’ 
dSphs at Qfint = 5 X ac- For this large integration angle, the cred¬ 
ibility intervals shrink when the halo size gets smaller, so that an 
underestimation of the halo size will also lead to an underestima¬ 



Figure Bl. Ratios of J-factors obtained with three different halo sizes Umax 
(in black, blue and red colours) to a reference J-factor with J?mj,x = 50 X 
r*. The smaller is Umax/Um®x' the more underestimated is the J-factor. 
The effect is stronger for core (dashed lines) than for cuspy DM profiles 
(continuous lines). 

tion of the credibility intervals. For the J-factors, this effect appears 
only at large integration angles (aint > olc). For the D-factors, the 
effect is more pronounced, and appears at all integration angles (not 
shown). 

Halo size estimation: comparison of the two methods In this 
work, we have considered two methods to estimate the halo size, 
detailed in Section 4.3. The first method uses the tidal radius as an 
estimation of the physical size of the halo; in the second approach, 
it is evaluated as the radius where the halo DM density and the 
MW DM density are equal. For a given dSph galaxy, the halo size 
is computed for each DM model accepted by the MCMC analysis. 
Figure B3 shows the distributions of Umax values obtained for Ursa 
Minor, using either the tidal radius estimation (red solid) or the 
equality of DM densities (dashed blue). Both distributions spread 
over more than one order of magnitude, with the mean of tidal radii 
distribution being systematically larger than the mean estimation 
from the DM density equality. This behaviour is found in all the 
dSphs. This trend is however not reflected on the astrophysical fac¬ 
tors, for which both methods give very similar results as shown in 
the bottom panel in figure B3. 


APPENDIX C: IMPACT OF CONTAMINATION FOR 
‘CLASSICAL’-LIKE DSPH GALAXIES 

Samples of stars from dSph galaxies may be contaminated by 
Milky Way and/or stream interlopers. Contamination in the con¬ 
text of Segue I will be discussed in Bonnivard, Maurin & Walker 
(in prep.). The less spectacular case of Fornax is presented below. 

Low impact contamination in Fornax Figure Cl shows the J- 
factors for Fornax from both the -weighted (red cross) and Pi- 
unweighted (cut Pi > 0.95, black circles) analyses. Contrarily to 
the other ‘classical’ dSphs, the two analyses give results which are 
in slight disagreement at the two-sigmas level (especially at large 
radii). The lower panels of figure Cl show the membership proba¬ 
bilities as a function of the projected radius R for Fornax (middle 
panel) and Carina (bottom panel): 7 out of the 8 ‘classical’ dSph 
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Figure B2. Distributions of J+95%CI/ jmed obtained for the 64 mock ‘ul- 
trafaint’ dSphs at aint = 5 x Qc- Three halo sizes are considered (short 
dashed blue, -Rmax = 5 x rj; long dashed green, Rmax = 20 x rj; solid 
black, Rm-dx = 100 x rj). At such large integration angles, smaller halo 
sizes reduce the credibility intervals. Underestimation of the halo size can 
then lead to underestimation of the CIs. 




Figure B3. Top: Distribution of halo sizes obtained for Ursa Minor, using 
either the tidal radius (solid red) or the DM densities equality (dashed blue) 
estimation. Both distributions show a large spread over one order of magni¬ 
tude, with the mean of the tidal radii distribution systematically larger than 
the mean estimation with the DM densities. Note that the mean light scale 
radius (r*) (obtained by fitting the surface brightness data) is typically one 
order of magnitude lower than halo sizes derived here. Bottom: Same as 
upper panel but for the J-factors. The systematic shift noted in the halo size 
does not propagate to the astrophysical factors. 


data display similar properties as Carina, i.e. most of the stars have 
Pi ^ 0.95 with velocities close to the mean value {v) for the ob¬ 
ject. Fornax shows a significant fraction of stars with intermediate 
Pi values (especially at large radii), whose velocity departs from 
the average. This is at the origin of the difference between the two 
reconstructed J-factors (top panel). Segue I shows a much stronger 
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Figure Cl. Top panel: J-factor as a function of the integration angle 
for the ‘classical’ dSph Fornax, for an analysis with and without weighting 
by membership probabilities (see text). The lower panels show the member¬ 
ship probability Pi of Fornax and Carina stars as a function of the projected 
radius R. The colour scale indicates the departure (red) from the mean ve¬ 
locity {v) of all stars (blue). 


dependence on the type of analysis, as will be discussed in Bon- 
nivard, Maurin & Walker (in prep.). 

Impact of contamination checked on mock data Jeans analy¬ 
ses on mock data with controlled levels of contaminant are helpful 
to investigate the impact and/or robustness of the reconstructed J- 
factors. Among the contaminant-free mock data described in Bon- 
nivard et al. (2015), we selected a couple of models (with rising, 
flat, or decreasing velocity dispersion, from a core and cusp DM 
profile) mixed with MW and stream stars (details will be presented 
in Bonnivard, Maurin & Walker in prep.). For each model, we cre¬ 
ated 1000 data sets with different levels of contamination from the 
MW and from the stream (denoted /mw and /stream respectively, 
with /mw + /stream ^ !)• Two Sample sizes were drawn to mimic 
‘classical’ (~ 300 — 3000 stars) and ‘ultrafaint’ (~ 30 — 100 stars) 
dSph galaxies. After reconshucting membership probabilities with 
the EM algorithm (Walker et al. 2009), we selected the contami¬ 
nated mock data that have at least 10% of their stars with interme¬ 
diate membership probabilities (0.1 ^ Pi ^ 0.95). The other sets 
of contaminated data correspond to a sample whose Pi are recon¬ 
structed with high certainty (similar to no contamination case). We 
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Figure C2. Analysis of reconstructed J-factors (Oint = Oc ~ 2r/i/d, 
large sample size) of contaminated mock data. The x-axis (resp. t;-axis) is 
the level of contamination from the Milky-Way (resp. stream) in the sample. 
The colour scale shows the ratio of the reconstructed to the true J value. The 
size of the symbols coirespond to the size of the Cl on the reconstructed J- 
factor. See text for the discussion. 


then ran the Jeans analysis, fixing the light and anisotropy profiles 
parameters to their true values to factor out non-essential ingredi¬ 
ents of the analysis. 

The results for the small sample sizes, in the context of 
Segue I’s analysis, will be presented in Bonnivard, Maurin & 
Walker (in prep.). Here, we present the result for the larger sam¬ 
ple sizes in relation to Fornax, and show that J-factors for such 
samples are much less affected by contamination. Figure C2 shows 
a comparison of the reconstructed to the true value for both the Pi- 
weighted and Pi-unweighted (cut Pi > 0.95) Jeans analyses on 
adversely contaminated mock data: the former analysis relies on 
the likelihood function equation (11), while the second takes care 
of removing all stars with Pi ^ 0.95. The maximal departure from 
the true value is observed close to maximum contamination (dashed 
line is /mw + /stream = 1). However, even in this worst case, the 
overshoot is within a factor 10 for most of the models. Four models 
show more important overshooting (up to a factor 10^), which are 
caused by a few misidentified contaminant stars with both large de¬ 
parture from the mean velocity (up to 100 km/s) and large member¬ 
ship probability. These stars could be easily identified and removed 
from a real data sample. We checked that none of the ‘classical’ 
dSph we studied had such extreme outliers. 

This leads to the conclusion that objects with a significant 
number of stars, as a result of having their membership probabil¬ 
ities robustly recovered (Walker et al. 2009), have a robustly re¬ 
constructed J-factor, even in the presence of large contamination. 
Regardless of the sample size, another conclusion is that whenever 
the fraction of stars with intermediate membership probabilities 
(0.1 Pi ^ 0.95) becomes significantly different from zero, the 
Pi — R—\Av\ plot is useful to identify likely contaminated objects: 
in the case of ‘ultrafaint’ dSph galaxies, the impact on the J-factor 
can be very important (Bonnivard, Maurin & Walker in prep.). 
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